Figure 5.5 Nusselt number distribution on the bottom plate of the channel 

with built-in circular tube 53 

Figure 5.6 Limiting streamlines on a plane close to the bottom plate of a 
channel with built-in circular tube and a delta winglet-pair in 
common-flow-up configuration. 54 

Figure 5.7 Streamline-plot on the horizontal midplane of a channel with 

built-in circular tube and delta winglet-pair in common-flow-up 
configuration 55 

Figure 5.8 Streamline-plot on the vertical midplane of a channel with 

built-in circular tube and delta winglet-pair in common-flow-up 


configuration 56 

Figure 5.9 Cross-stream velocity vectors at different x-locations behind the 

delta winglet-pair 57 

Figure 5.10 Cross-stream velocity vectors at different x-locations behind the 

rectangular winglet-pair 58 

Figure 5.11 Streamlines on the cross-stream planes at different x-locations 

behind the delta winglet-pair (b/D=1.67) 59 

Figure 5.12 Streamlines on the cross-stream planes at different x-locations 

behind the rectangular winglet-pair (b/D=l .67) 60 

Figure 5.13 Streamlines on the cross-stream planes at different x-locations 


behind the delta winglet-pair (b/D=1.08) 


61 



Figure 5.14 Three dimensional particle paths behind the delta winglet pair 

in common-flow-up configuration (b/D= 1 .67) 62 

Figure 5.15 Three dimensional particle paths behind the rectangular winglet-pair 
in common-flow-up configuration (b/D= 1 .67) 63 

Figure 5.16 Three dimensional particle paths behind the delta winglet-pair in 

common-flow-up configuration (b/D=1.08) 64 

Figure 5.17 Nusselt number distribution on the bottom plate of a channel with 
built-in circular tube and delta winglet-pair in common-flow-up 
configuration 65 

Figure 5.18 Nu s for a plane channel, channel with built-in circular-tube and 


channel with built-in circular-tube and delta winglet pair 66 

Figure 5.19 Span-averaged Nusselt number distribution for different Reynolds 

numbers 67 

Figure 5.20 Nu s for two different types of winglet-pair 68 

Figure 5.21 Nu s for the delta winglet-pairs with different b/D 69 

Figure 5.22 Grid independence test 70 

Figure 5.23 Comparison between the present computation and experimental 


results 


71 



N omenclature 


B channel width 

b lateral distance between trailing edges of a winglet-pair 
C p specific heat of the fluid 

CV control volume 

D diameter of the cylinder 

F mass flux through a cell face 

H channel height 

k thermal conductivity of the fluid 


Nu 

Nu , 


local Nusselt number, 


r d0 

< 


1 w 


(l-6?J \dZ 

span averaged Nusselt number, J Nu dy / jdy 


p static pressure 

Re Reynolds number based on channel height ^ m — 

M 

S surface area of a cell face 
t time 

u axial velocity 

v spanwise velocity 

w normal or vertical velocity 

x axial dimension of coordinates 

y spanwise dimension of coordinates 

z normal or vertical dimension of coordinates 

T temperature 


Greek Letters 

a i, a 2, a 3 coefficients of linearly independent unit vectors 



p 

e 


density of the fluid 

non-dimensional temperature (T -T„) /(T w - ) 


Subscripts 

b bulk condition 

j cell face 

i cell center 

w wall 

Superscripts 

n, n + 1 time level 



Chapter 1 


Introduction 


1.1 Motivation 

Fin-tubes are commonly used in gas-liquid crossflow heat exchangers. In such 
heat exchangers, the gas generally flows across the tubes and the liquid flows 
inside the tubes. The fins act as extended surfaces providing the bulk of the heat 
transfer area on the gas side. Even with the extended surfaces, the dominant 
thermal resistance is on the gas side. In order to achieve significant enhancement 
of heat transfer on the gas side, strategies must be developed. The heat transfer 
coefficients on the fin and the outer surfaces of tubes are to be enhanced keeping 
the pressure penalty within a modest limit. 

1.2 The Importance in Technology and Applications 

Heat exchangers are widely used in power and process industries. Improvement in 
their performance is of great technical, economical, and ecological importance. 
Existing air-cooled condensers in geothermal power plant use fin-tube heat 
exchangers with circular tube. The size of heat exchangers is quite large, 
consuming a large fraction of overall capital cost of these plants. Accomplishment 
of more efficient heat transfer on the airside of the air-cooled condensers can 
reduce its size and can improve the overall performance of all the geothermal 
plants. A modest additional pressure drop across the heat exchangers must support 
the enhancement of heat transfer. Figure 1.1 shows a schematic diagram of the 
core region of fin- tube heat exchangers. The protrusions on the flat surfaces can 
disrupt the growth of thermal boundary layer and increase the transport 
coefficients. The protrusions in the form of delta wings or winglets can bring 
about the desired augmentation at the expense of small increase in pressure drop. 
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1.3 Focused Area and Specific Outline of the Present Problem 

One of the effective methods of reducing the thermal resistance is induction of 
longitudinal streamwise vortices in the flow field. The longitudinal vortices 
interact with an otherwise two-dimensional boundary layer and produce a three- 
dimensional swirling flow that mixes the near wall fluid with the free stream. 
Such kinematic mixing strongly enhances the entrainment of fluid from the 
periphery to the core region of the flow field. As a result, thermal boundary layer 
is disrupted and heat transfer is enhanced. The longitudinal vortices can be 
created by mounting delta wing or delta winglet type of vortex generators on the 
flat surface (Figure 1.2). The longitudinal vortices develop along the side-edges of 
the delta-winglets due to the pressure difference between the front surface facing 
the flow and the back surface. These vortices are also called the streamwise 
vortices since the axes of vortices are aligned to the flow direction. The pressure 
drop caused due to such vortices is modest because form drag for such winglet- 
type slender bodies is low. Vortex generators can be mounted behind the tubes in 
fin-tube heat exchangers in two common configurations - (i) common-flow-down 
and (ii) common-flow-up as proposed by Pauley and Eaton (1988). In a common- 
flow-down configuration, the transverse distance between the leading edges of the 
winglet pair is less than the transverse distance between the trailing edges. On the 
other hand, in a common-flow-up configuration, the transverse distance between 
the leading edges of the winglet pair is more that between the trailing edges. 

In order to achieve the objectives enumerated earlier, present investigation has 
been carried out for the flow and heat transfer in fin-tube cross-flow heat 
exchangers with delta-winglets mounted behind the tube in common-flow-up 
arrangement. Numerical investigation on the related topic has been carried out 
earlier by Biswas et al. (1994). They have used the winglets with common-flow- 
down arrangement. Torii et al. (2001) have proposed a novel technique of placing 
the winglets with common-flow-up arrangement. Figure 1.3 presents a sectional 
view of the suggested arrangement for the fin-tube heat exchangers with delta- 
winglets placed behind the tube in common-flow-up configuration. For such a 
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configuration, the full Navier-Stokes equations together with the energy equations 
are solved and a detailed analysis of flow structure together with heat transfer 
characteristics is carried out 

1.4 Layout of the Thesis 

In Chapter- 1 we have already introduced the genesis of the problem. Chapter-2 of 
the thesis provides a review of the literature relevant to the understanding of the 
basic mechanism involved in augmentation of heat transfer and the algorithms 
related to solving the Navier-Stokes equations. The mathematical formulation of 
the problem is presented in Chapter-3. Here, the geometry of interest, the 
governing equations, boundary conditions are discussed. Chapter-4 discusses the 
grid generation technique, the solution algorithm and the mathematical foundation 
of the finite volume formulation in detail. Chapter-5 discusses the results of the 
present investigation with emphasis on flow physics and associated heat transfer. 
Chapter-6 includes the concluding remarks and the scope for the future research. 
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Figure 1.1 Schematic diagram of core region of a fin-tube heat exchanger 
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Figure 1 .2 Delta-Winglet type vortex generators on flat surface 








Chapter 2 


Literature Survey 


2.1 Introduction 

This chapter is divided in three different sections. The first section provides an 
overview of work done by various researchers in the area of enhancement of heat 
transfer. Quite appreciable amount of experimental and computational research 
have been carried out in the area of heat transfer augmentation. Especially 
concerned literature on enhancement through manipulation of surface geometry, 
which has been analyzed by many researchers, is presented in this chapter. The 
next section focuses on the work done m the field of incompressible flow 
modeling and associated heat transfer problems. The third section discusses 
aspects of unsteady flow and turbulence. Finally the benefit of current problem is 
established from the existing literature. 

2.2 Augmentation of Heat transfer 

The present aim of the study is to enhance the heat transfer on the gas side of the 
fin-tube heat exchangers and to minimize the space occupied by the equipment for 
a desired rate of heat transfer. With this intent we would like to study different 
investigations related to the augmentation of heat transfer suitable for the 
applications mentioned above. 

The representative element of a fin-tube heat exchanger consists of a rectangular 
channel with built-in circular tube (Figure 2.1). In channel flows, rate of heat 
transfer between fluid and channel walls deteriorates as boundary layer grows on 
the channel walls and flow tends to become fully developed. Protrusions can be 
mounted on this channel walls in order to disrupts the growth of boundary layers 
and thereby enhance the heat transfer between the flowing fluid and channel wall. 
To achieve this goal, protrusions in form of slender delta wings or winglets can be 
mounted . The longitudinal vortices are generated along the side edge of the wing 
shaped vortex generators due to pressure difference between the front surface 
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facing the flow and the back surface. Investigations on the enhancement of heat 
transfer by vortex generators have been carried out in the last two decades by 
Fiebig et al. (1986), Turk and Junklan (1986), Eibeck and Eaton (1987), Fiebig et 
al. (1989), Tigglebeck (1994), Yanagihara and Toni (1990), Fiebig et al. (1990), 
Fiebig et al (1991), Sanchez (1989), Biswas et al (1994) and Valencia (1996). 
The idea of using vortex generators in compact heat exchangers has already been 
implemented by industry. Fiebig and his co-workers (Fiebig et al. 1986) 
experimentally observed enhancement in heat transfer using vortex generators on 
the primary heat transfer surfaces. Subsequently a host of other investigations 
were taken up by different groups and the use of vortex generators evolved as a 
device for enhancement in heat transfer. Biswas and Chattopadhyay (1992) have 
considered the influence of Reyonlds number and angle of attack on skin friction 
and Nusselt number for delta wing type vortex generators in a channel. 

Fiebig and his co-workers (1991) conducted an experimental study to compare the 
performance of delta wings, delta winglets, rectangular wings and rectangular 
winglets as the enhancement device. Delta wings and winglets were found to be 
more effective than rectangular wing and winglets. For the simple case of a flat 
plate with a delta wing attached to its leading edge, Gentry and Jacobi (1997) 
presented a method for selecting the preferred wing aspect ratio and attack angle. 
The interaction of vortices with the boundary layer and its effect on heat transfer 
is a subject of interest to many researchers. Some notable investigations include 
Taylor-Gortler vortices in boundary layers on concave surfaces, the horseshoe 
vortices formed by an obstruction protruding from the inner surface and wingtip 
vortices impinging on a downstream surface. The embedded vortex is capable of 
strongly perturbing the boundary layer and influencing the heat transfer 
characteristics. In addition, longitudinal vortices usually maintain their coherence 
over a long streamwise distance. As a consequence, the heat transfer effects 
behind the vortex generators are very persistent. Eibeck and Eaton (1987), 
Westpal and Mehta (1987) have worked extensively in this field with more focus 
to the turbulent boundary layers. Eibeck and Eaton (1987) have conducted 
experiments on longitudinal vortices embedded in a turbulent boundary layer and 
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resultant heat transfer effects. Longitudinal vortices are found to influence the 
heat transfer behavior significantly. Local Stanton number increase is as high as 
24 percent resulting in a net increase in spanwise average heat transfer coefficient. 

2.3 Literature on Incompressible Flow Modeling and Associated Heat Transfer 
Problems 

During past three decades various numerical techniques for modeling 
incompressible flows have been developed. Such techniques involve numerical 
solution of full Navier- Stokes equations and energy equations for obtaining the 
temperature field. The major difficulty encountered during the solution of 
incompressible flows arises from the absence of any explicit equation for pressure 
and due to the nature of spatial coupling of the pressure and velocity. For 
incompressible flow problems, pressure does not have its usual 
thermodynamically meaning. Here it is a relative variable, which adjusts itself 
instantaneously for the condition of zero mass divergence to be satisfied at all 
computational cells. This behavior is reflected from the fact that in an 
incompressible fluid speed of sound is infinite. As a consequence, the pressure 
field cannot be solved by explicit time advancement procedure, instead it requires 
at least a partially implicit determination which takes into account the coupling 
between pressure and velocity fields as well as the effects of velocity boundary 
conditions. This aspect is the most distinctive feature of the primitive variable 
formulation of the incompressible flow modeling. The difficulty of determination 
of pressure field can be resolved in the stream function-vorticity approach. But 
this approach loses its attractiveness when three-dimensional flow is computed 
because of absence of a single scalar stream function in the three-dimensional 
space. Efforts has been made so that two as well as three-dimensional problems 
could be computed following a primitive variable approach without encountering 
non-physical wiggles in pressure distribution. As a remedy, it has been suggested 
to employ a different grid for each variable. Harlow and Welch (1965) have used 
a staggered grid for the dependent variables in their well-known MAC (Marker 
and Cell) method. The MAC method of Harlow and Welch is one of the earliest 
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and widely used explicit methods for solving the full Navier-Stokes equations. 
The MAC method uses a layer of imaginary cells around the boundary of the 
physical domain necessitating the updating of the boundary conditions after every 
change in the internal velocity and pressure values. The onginal version of MAC 
method has been modified by Harlow and Amsden, popularly known as 
Simplified MAC (SMAC, 1970). Nicholas and Hirt (1971) and Hirt and Cook 
(1972) have also used a modified version of MAC algorithm for application of 
free surface flows. The MAC method has been successfully used by many 
researchers to simulate even highly unsteady and turbulent flows (Robichaux, 
Tafti and Vanka, 1992). It has been experienced that the MAC method is indeed 
very efficient in the study of temporal flow development. It has stability 
restrictions on the time increment, which slows down the calculation for steady 
flow considerably. Since implicit methods have no such restrictions, they are 
more attractive. 

Mukhopadhayay et al. (1993) have developed a numerical method for predicting 
viscous flows for incompressible geometries. Integral mass and momentum 
conservation equations are deployed and these are discretized into algebraic form 
through numerical quadrature. The physical domain is divided into a number of 
non-orthogonal control volumes that are isoparametrically mapped onto standard 
rectangular cells. Verma and Eswaran (1996) and Verma and Eswaran (1997) 
have developed a finite volume based flow solver for complex geometries using 
overlapping control volumes. 

Patankar and Spalding (1972) have introduced an efficient method known as 
SIMPLE (Semi-Implicit Method for Pressure Linked Equations). This method is 
based on finite volume discretization of the governing equations on a staggered 
grid. In order to improve the convergence involved in the solution of pressure 
equation, several variants of SIMPLE algorithm have been developed. The 
SIMPLER algorithm of Patanker (1981) and SIMPLEC algorithm of Van 
Doormaal and Raithby (1984) are improvements on SIMPLE. Although changes 
to incorporate SIMPLEC into SIMPLE algorithm are minor, the consequences 
can be great as it eliminates the approximations made in SIMPLE while deriving 
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the equation for pressure-velocity correction. Garg and Maji (1987) have applied 
SIMPLEC (Van Doormaal and Raithby, 1984) method for solution of viscous 
flows through periodically converging-diverging tubes. Application of finite 
volume methods using non-orthogonal coordinates and collocated grid 
arrangement is reported by Rhie and Chow (1983) and Peric (1985). Peric et al. 
(1988) have shown that collocated arrangement converges faster than the 
staggered variable arrangement and has advantages when extensions such as 
multigrid techniques and non-orthogonal grids are considered. An appropriate 
pressure interpolation avoids the pressure splits in collocated grid arrangement. A 
finite volume based procedure with momentum interpolation has been 
successfully developed by Majumdar (1988) using such collocated velocities and 
pressure. Kobayashi and Pereira (1991) have modified the momentum 
interpolation method suggested by Peric (1985) and named it as Pressure- 
Weighted Interpolation Method Corrected (PWIMC). In this method the non- 
orthogonal terms in the momentum equations were solved explicitly, whereas in 
the pressure corrections they were dropped. 

It is clear that substantial progress has been made for the development of 
algorithms for complex incompressible flow simulations but none of the 
prescriptions is universal. Depending upon the nature of flow and geometry etc., 
one can always go for the best-suited discretization procedure for solving time- 
dependent three-dimensional Navier-Stokes equations. 

2.4 Unsteady flows and Turbulence 

In the experimental investigation of transport process, Cantewell and Coles 
(1983) have studied the near- wake of a circular cylinder at a very high Reynolds 
number. The periodic and random components have been shown to have 
comparable amplitudes. A considerable emphasis has been placed on the topology 
of the unsteady mean flow. This emerges as a pattern of centers and saddles in a 
frame of reference moving with the eddies. The kinematics of the vortex 
formation process has been described through the critical point theory. The 
important conclusion of the work is that the turbulence production is concentrated 
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near the saddles. Entrainment of the fluid into the wake is also found to be closely 
associated with the formation and location of the saddles. 

The wake behind circular tube shows unsteady characteristics even at low 
Reynolds numbers. Reynolds number behaves as control parameter in a 
dynamical system and as it is increased, detachment of the free shear layer and the 
consequent shedding of vortices, generation of harmonics, three-dimensionality in 
nominally two-dimensional geometries, onset of chaos and transition to 
turbulence results. This by itself is a sequence of transitions, each being identified 
with a critical Reynolds number and an associated flow topology. A non-linear 
dynamical system can go to a chaotic state through three distinct routes. These 
involve period doubling of frequencies (Feigenbaum, 1980), the Ruelle-Takens- 
Newhouse route through quasi-periodicity (Ruelle and Takens, 1971) and 
intermittency (Manneville and Pomeau, 1980). Vittori and Blondeaux (1993) have 
reported a quasi-periodic route to chaos in their numerical study of two- 
dimensional oscillatory flow around a circular tube. They found that the system 
gets phase-locked before it becomes chaotic. In the quasi-periodic route, the 
dynamical system initially at a steady state, becomes unsteady as the control 
parameter crosses a particular threshold limit. The nonlinearity of the system 
generally increases as a result of an increase in control parameter. For example, 
an increase in the Reynolds number diminishes diffusion and dissipation 
mechanisms, resulting in greater prominence of the non-linear acceleration terms. 
As a consequence of instability, the dynamical behavior approaches a limit cycle 
caused by a Hopf bifurcation. Similarly, further increase of Reynolds number 
leads to appearance of three frequencies. Guzman and Amon (1996) have shown a 
quasi-periodic and frequency locking route to chaos in two as well as three- 
dimensional converging-divergence channels. They have found that the flow 
undergoes a rapid transition to a chaotic three-dimensional state at a Reyonlds 
number of 500 starting from a laminar two-dimensional state at a Reyonlds 
number of 200. In the work of Pulliam and Vastano (1993), the transition to chaos 
of an open unforced two-dimensional flow past an aerofoil has been documented. 
The salient findings of their study are that the system undergoes period-doubling 
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bifurcations to chaos as the Reynolds number is increased from 800 to 1600, The 
dominant harmonics present in the signal can be identified with the help of FFT 
(Fast Fourier Transformation) of the time series data. 

2.5 How the Current Problem Benefits from Existing Literature 

The numerical investigations of Biswas et al. (1994) and the experimental 
findings of Valencia et al. (1996) reveal the effective utilization of vortex 
generators for the purpose of enhancement of heat transfer while producing less of 
a pressure drop. In the above investigations, the enhancement of heat transfer 
from the fin surfaces is achieved by placing delta-winglet type vortex generators 
ion the flat fin surfaces in the neighborhood of the tube. In the related study of 
Tsai et al. (1999), the use of wavy-type extended fin surfaces as device for 
enhancement of heat transfer in the fin-tube heat exchangers has been deployed. 
The longitudinal vortices or streamwise vortices get developed along the side 
edges of the delta-winglets due to pressure difference between the front surface 
and the back surface. These vortices interact with the boundary layer and produce 
a three-dimensional swirling flow that mixes near wall fluid with the free-stream. 
The mechanism strongly enhances the entrainment of fluid from the periphery to 
the core region of the flow filed. Thus the thermal boundary layer is disrupted and 
the local heat transfer rate is enhanced. The additional pressure losses are modest 
because the form drag for such winglet-type slender bodies is low. 

The air-cooled condensers of geothermal plants consist of the same fin-tube 
arrangement as shown in Figure 2.1. Air is forced through several rows of these 
fin-tubes by large fans. The condensers units can be very large, consuming a large 
fraction of the overall capital cost of these plants. In the present research, a three- 
dimensional numerical model has been formulated to provide a better insight in to 
the flow physics. The complete Navier-Stokes equations together with governing 
equations of energy are solved in an element of heat exchanger as shown in 
Figure 1.3. The heat exchanger element consists of a rectangular channel with a 
built-in circular tube and a winglet pair. A detailed analysis of the flow structure 


12 



along with heat transfer characteristics in such an element is studied in this 
investigation. In general, the delta-winglet type vortex generators can be mounted 
on the flat surfaces with common-flow-down configuration. For such a 
configuration, the transverse distance between the leading edges of the winglet 
pair is less than the transverse distance between the trailing edges of the winglet 
pair (Pauley and Eaton, 1988). In the present case of common-flow-up 
configuration, transverse distance between the leading edges of the winglet pair is 
more than that between the trailing edges. Torii et al. (2000) have indicated that 
the delta winglet pair with common-flow-up configuration creates constricted 
passages in the aft region of the tube, which brings about separation delay. The 
fluid is accelerated in the constricted passages, and as a consequence, the point of 
separation travels downstream. Narrowing of the wake and suppression of vortex 
shedding are the obvious outcome of such configurations. Since the fluid is 
accelerated in this passage, the zone of poor heat transfer is also removed. The 
poor heat transfer zone is created in the wake of the tube in absence of any vortex 
generators. The common-flow-up configuration of the delta winglet pair emerges 
as the most attractive enhancement technique for the air-cooled condensers. 
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Figure 2.1 Element of existing fin-tube heat exchanger 


Chapter 3 


Mathematical Formulation 


3.1 Introduction 

Arising out of practical importance of enhancement of heat transfer in fin-tube 
heat exchangers a need is felt to analyze the flow field in a representative element 
of the heat exchanger. The computation involves determination of flow and heat 
transfer characteristics in a horizontal channel with a built-in circular tube and a 
pair of winglets. The numerical simulation can provide point by point information 
of flow structure and temperature distribution in an element of a fin-tube heat 
exchanger along with the effect of various parameters on heat transfer 
enhancement. 

3.2 Statement of the Problem 

Figure 1.3 represents the computational domain. Two neighboring fins form a 
channel of height H, width B = 11.25H and length L = 20H. The circular tube of 
diameter D = 4.0H is located at a distance Li= 1.59D from the inlet. The blockage 
ratio, D/B has been kept 0.36. The position of the winglet is given by X A = 7.6H, 
X b =9.8H, Y a = 1.7H, and Y B =3.66H. The height of the winglet, h=0.347H and 

the angle of attack, (3=138.61 ° Air has been used as the working fluid; hence the 
Prandtl number of this study is 0.7. 

The flow field in the channel with built-in circular tube and a winglet pair is 
characterized by the following parameters 

• Reynolds Number 

• Position of the oval tube in the channel (Li/L) 

• Velocity profile at the channel inlet 

• Channel height (H/B) 
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• Angle of attack of the winglets ((3) 

• Height of the winglets (h/H) 

® Position of the winglets (X A /L, Y A /L) 

3.3 Governing Equations 

The three-dimensional Navier-Stokes equations for laminar flow of an arbitrary 
spatial control volume V bounded by a closed surface S can be expressed in the 
following general convection-diffusion-source integral form 


3 _ 

dt 



+ 


J/?u-dS 


s 


= 0 


(3.1) 



+ jLou^-r,vrf-ds =\s,d.v 

S V 


(3.2) 


Where p represents the fluid density, u is the fluid velocity, <f> stands for any 

vector component or scalar quantity, S 0 is the volumetric source term. For 

incompressible flow of Newtonian fluid, the equation takes the form 

Ju.d,S=0 (3.3) 

£ 

|-J>0dV + l[pu<f>-r^<p].dS = js,dV (3.4) 

Ot V s v 

and the source term for momentum equation becomes — — f p I.dS where I is the 

P J s 

unit tensor. In this formulation we work with Cartesian components of velocity. 
So 0 can be the three-cartesian components of velocity u, v, w as well as any 
scalar e.g., temperature which need to be determined. The variables of the general 
transport equation are given in Table-3.1. 
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Table 3.1: Values of the variables in the general transport Equation (3.2) 


Equation 

0 

r, 

S<t> 

Continuity 

1*" 

1 

0 

0 

Momentum 

Uj 

F 

dp 

dxj 

Energy 

T 

k 

7Tp 

0 


3.4 Boundary Conditions 

The governing differential equations are elliptic in space and parabolic in time. 
We need the boundary conditions for all the confining surfaces. The boundary 
conditions of interest in the present investigation are: 


Top and bottom plates 

u = v = w = 0 , (No Slip boundary condition) and — — = 0 

az 


T = T W ( T w represents wall temperature) 


Side Wall 


d u _ 3 w 
By By 


v = 0 


dT 

dy 


= 0 


(Free Slip boundary condition) and 


d P 

dy 


= 0 


Channel Inlet 
u — , v=w=0, 


and 


3 £ 

d x 


= 0 


T = T„ 

Channel Exit 
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The mass flux through the outlet boundaries is found by means of a continuative 
outflow condition (Orlanski (1976)), which allows changes inside the flow field 
to be transmitted outward, but not vice-versa. 

d <p d (j) 

-r— + ” m v - = 0 (Where <p represents u, v, w or T) 
at ox 

P = P~ 

Obstacles: Circular tube and the Winglet pair 

r) n 

u = v=w= 0, and — — = 0 (where n signifies normal direction) 
an 

and T = T w 
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Chapter 4 


Grid Generation And Method Of Solution 


4.1 Introduction 

In order to solve numerically the governing equations in the prescribed 
computational domain with given initial and boundary conditions, two major 
discretization steps have to be performed. First being the discretization of domain, 
numerical grid is generated for computation Second is the discretization of the 
governing equations. Various methods of discretization can be used for a given 
set of partial differential equations. The most common methods used are finite 
differences, finite volumes and finite elements. In present computation, a finite 
volume scheme has been used for spatial discretization of governing equations. 
Conservative formulation of discretization is followed, assuring automatically the 
satisfaction of the global balance of the conserved quantities, independent of the 
coarseness of the numerical grid. 

4.2 Grid Generation 

Figure 4.1 shows the schematic three-dimensional grid for the computation. A 
differential equation is used to generate the body fitting type two-dimensional 
grid where circular tube and a pair of winglets are placed. The two dimensional 
grid is stacked layer by layer in third direction to generate the grid in that 
direction. 

Initially two-dimensional grid is generated by the transfinite interpolation, which 
uses essentially linear interpolation scheme to compute the interior points by 
using the points from the computational domain boundaries. Thus an algebraic 
grid is obtained, which is further smoothened by the use of Partial Differential 
Equations (PDE). In this technique, a system of elliptical PDEs is solved for 
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obtaining the location of grid points in the physical space, whereas on the 
computational space, transformed grid is of rectangular shape with uniform 
spacing. A Laplace equation or a Poisson equation with a Dirchlet boundary 
condition can be used for this purpose. Standard Poisson equation is of the form 

V 2 £ = P 

V 2 ti = Q (4.1) 

Where P and Q are known as control functions. Here we use Laplace equation 
in the following way 

V 2 £ = 0 

V 2 r) = 0 (4.2) 

The above equations can be solved by finite difference technique to get the 
location of the interior points from the boundary values. The algebraic mapping 
can generate initial guess, needed to solve the Laplace equation. The Laplacian 
operator can provide quite smooth grids. The grid lines generated remain equally 
spaced if boundary curvature is absent but show a tendency to concentrate near 
curved boundaries. By using the appropriate control functions in the Poisson 
Equation, this drawback can be removed. But the use of control functions add 
more complexity in the transformed domain. 

4.3 Computational Procedure for Grid Generation 

The algorithm used for the generation of the present grid can be outlined as 
follows: 

1 . Geometrical data is used as input. 

2. The computational grid is defined on £,-r| plane based on the number of grid 
points. 

3. Boundary points along the edges of the physical domain are defined. 

4. Algebraic grid is used to generate the initial guess by using liner interpolation 
technique, which is used by the Laplacian operator. 

5. Using four-point formula, guess is further improved a little. 
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6. Elliptic grid generation system constitutes of following equations derived 
from the Laplace operator 

(X " 2y(? -V^r) ~ ^ ' (4.3) 

ocyii, - + Xy m = 0 

Where 

j 

CC = x n 2 + y n 2 
(3=^x n + y^ y n 

W+n 2 

Here subscripts indicate partial derivatives of Cartesian variables with respect 
to curvilinear variables. 

These coefficients are determined using finite difference approximations. Their 
values in these expressions are provided by the initial distributions for the first 
iteration, and subsequently from the previous iteration, i.e., the computation of 
coefficients lag by one iteration level. The iterative solution continues until a 
specified convergence criterion is met. 

Because of discretized approximation of differential equations, the truncation 
error in numerical analysis is much dependent on the nature of gird. We target 
for a grid, which gives minimum truncation error. 

The Parameters that affect the quality of the grid are: 

1. Transformation Jacobian 

2. Skewness 

3. Aspect Ratio 

4 . Adjacent cell ratio 

A desirable grid should have nonzero transformation Jacobian almost at each 
location of the grid cells. It should be of similar order of magnitude with no quite 
small values appearing anywhere. It should possess least possible skewness, 
aspect ratio close to 6 and adjacent cell ratio close to one. 
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4.4 Finite Volume Method 


The conservation equations are discretized employing the finite-volume approach 
of Eswaran and Prakash (1998). The solution domain is divided into a number of 
contiguous (finite) control volumes (CV). The coordinates of the control volume 
vertices are calculated by the grid generation procedure. These vertices, which are 
assumed to be connected by straight lines, define the control volume. A 
collocated grid arrangement is employed and all the dependent variables (u, v, w, 
p, T) are defined at the same location, the centroid of the control volume 
(Figure 4.2). Cartesian velocity components have been used. The symbols, E, W, 
N, S, T and B indicate the six neighboring control volume centers for the east, 
west, north, south, top and bottom neighbors. The face center points e, w, n, s, t 
and b are located at the intersection of the lines joining the midpoints of the 
opposite edges. Edge centers te, be, ne and se are the midpoints of the edges that 
form the east face with e as the center of the cell face. Similarly tw, bw, nw, sw, 
tn, bn, ts & bs are located at the corresponding center of edges that form west, 
north & south face. 

The advantages of the collocated grid over staggered grids, especially when the 
non-orthogonal coordinates are: 

• There is only one set of control volume as all the variables share the same 
location. 

• The convection contribution to the coefficients in the discritized equations is same 
for all variables. 

• When coordinate oriented velocity components are used, it results in simpler 
equations. 

• Number of constraints on the numerical grid is reduced. 
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4.5 Surfaces areas and Volumes 


Figure 4.2 shows the comer points of the control volume numbered as 1,2,3, etc. 
The outward surface normal and volume can be found in the following manner as 
suggested by Kordulla and Vinokur (1983). 

Defining r tJ = iy - where r, and rj indicate the position vectors of points i and j 
respectively, we have the formulate the expressions for calculating the surface 
area and volume for different faces. 


East Face 

S e ~ ~ (r 7 4 x r 83 ) 

(4.4) 

West Face 

S w = ~ (rig x r 52 ) 

(4.5) 

North Face 

S n = ^ (r 27 x r 63 ) 

(4.6) 

South Face 

S s = i (ris x r 4S ) 

(4.7) 

Top face 

S t = j (r 75 x r 68 ) 

(4.8) 

Bottom Face 

Sb= ^(ri 3 xr 24 ) 

(4.9) 


The volume of the cell is calculated from the cell coordinates, with the 
assumption that the linear segments to form the six cell faces join the cell corners. 
The volume is calculated as follows 

V = |r 71 .(S s + S b + S w ) (4.10) 


4.6 Discretization Procedure 

The main steps for the discretization of the integral governing equations involve 
the calculation of convection and the diffusion fluxes and source terms. The rates 
of change and source terms are integrated over the cell volume, whereas the 


23 



convection and the diffusion terms form the sum of fluxes through the faces of the 
control volume. 


4.6.1 Discretization of Continuity Equation 

Equation 3.1 is discretized in the following manner: 

l yOU-dS = 5>(u-S), = 2>u,-S, (4.11) 

j=e,w,n,s,t,b j 

Where Sj is the surface vector representing the area of the j th cell face and uj is the 
velocity defined at the face center j. 

In discretized form continuity equation becomes 

SF J =F e +F w + F„ + F, + F,+fi,=0 (4.12) 

J 

where the Fj is the outward mass-flux through face j, defined by: 
p j = P • Sj 

4.6.2 Discretization of General Equation 
(a) Rate of Change 


The value of dependent variable (p at the centroid of the control volume represents 
an average over the CV as a whole. Thus 


i-r P ^v M v ) ? r{pE )i. apV t J'-r, 

dt^ v A t A t 

Where V is the volume of the cell. 


(4.13) 
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(b) Convection Fluxes 


The surface integral over convection flux of variable 0 can be approximated in 
the following form 


J s Pu^dS - 5>0>-S), =E F ^, 

J J 


(4.14) 


Where (f) ] is the value of q> at the center of face j. Thus 


j s P*<P • dS= F e 0 e + FJ W + F n 0„ + F S 0 S + F t 0, + F b <p b (4. 15) 

Where, 0 e is the (interpolated) value of the variable 0 at the east face center, etc. 
This can be evaluated by using a central difference linear interpolation between 
the neighboring nodal values 0 P md0 E . At east face the value of (j. ) e is given by 


V, 




0 = — — — d)p + — — — 0 r 

v E +v P ^ p v e +v/ e 


(4.16) 


Where V E and V p are volumes of the cells around the points E and P respectively 
and 0 E and0 p are the values of the dependent variables at these points. In a 
collocated grid system, all dependent variables are defined at the same location 
hence exactly the same interpolation scheme is used to express all of them at the 
interfaces. The central difference approximation to compute the convection flux 
may lead to numerical stability therefore it is blended with fist order upwind 
differencing scheme (UDS) as follows 

F,&=( F,h ) CDS +Y [ ' ( W, f DS 1 < 4 - 17 ) 

The upwind differencing scheme is based on assumption that the convected flux 
at cell face is equal to that at upstream cell along the same coordinate direction. 
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Thus the value (p e at east face is assigned the value if u e > O.This can be 
conveniently summarized as 


Fe <t> e -F e ( 77 - <p P + 77 — v, <P E )+y{ (j>p [ [F e ,o] ]-^e[ [-F f ,0] ]- 


V'+V, 


FVp 


\ K 


Fe(~hr<P P + 77 ^^) } 


Vn+V, 


v E +v„ 


(4.18) 


The symbol [ [ , ] ] signifies greatest of the two quantities inside the brackets and 
y is upwind factor. Similar expressions can be written for the other cell faces 
(Eswaran et al. 1998). 

-4)+t{ trl [f«0] ]-#„[ [-F„0] ]- 


F 1 v<pw~F t v( 


Vu, . y o 


— <P P + — 

v w +v P ' V ut +V B 


Y W n Y P 


V, v. 


+ 77 ^^) } 

Y W + V P 


v w +v p 


(4.19) 


F„^F„(-^— A. + ^)+y{ M [F,,0] ]-0 N [ [-F,,0] ]- 

V N -rVp V N + V ^ 


F„ ( — <p P + — £ ^ ) } 

y*,+y» y„+K 


(4.20) 




' N ' * P 


V. 


v D 


Fs<f>s=F s ( 77 ^* 7r </> p + 77 ^ 7r <f> s ) +y{ <f>p [ [Fj,0] ]-0s[ [- F it 0] ]- 


y s +y P r " y s +y P 


v. 


V, 


F s ( — <p P + — £- 

Vc+v» K+y 


A) 1 


(4.21) 


F,<^F,(-^— (f> P + 0 r ) ^p[[F,0]]-^ t [[-F,0]]- 

r TT ' r r 


y. . y D 

— 'Pp 


F 1 ( —- J —— <Pp + ) } 


yr+y P rr y r +y P 


(4.22) 


Fb</>b= F b ( 


y„ 


^[[^O]]-0 B [[-F fc O]]- 

*bFv p 


y B 


F|,( v,+v, + v b +v/ b) * 


(4.23) 
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In a fully implicit method the upwind parts of the above equations are implicit and 
they are incorporated in the coefficients of the unknown velocity during the 
pressure-velocity iterations. The CDS terms on the other hand are evaluated using 
the previous iteration values and used as a source term on the right side of the 
same equation. This is the so-called deferred correction approach of Khosla and 
Rubin (1974). Multiplication of the explicit part by a factor of y (0 < y < 1) 
allows the introduction of numerical diffusion. For the first order upwind 
difference scheme, y = 1 and for the second order central difference, y = 0. The 
deferred correction approach enhances the diagonal dominance of the coefficient 
matrix. However, the present solution scheme being explicit, the accuracy 
depends on the value of y. 

(c) Diffusion Fluxes 

Peric(1985) calculated the expression for the diffusion fluxes in the finite volume 
approach. However, the exact expressions in three dimensions are not explicitly 
available in open literature. Eswaran and Prakash (1998) have evaluated the 
diffusion flux of a variable <j> through the cell faces in the following manner 

| s r,V«>-dS= 2 (4.24) 

j=e,w,n,s,t,b j 

For any face, Surface vector can be expressed as 

S j = cqn 1 + a 2 n 2 + a 3 n 3 (4.25) 

Where n\ n 2 and n 3 are any three linearly independent (not necessarily 
orthogonal) unit vectors. Therefore 

V0-S. = V0- (a,n‘ + o,n 2 +tf 3 n 3 ) (4.26) 

= a { V^-n 1 + a 2 V0-n 2 +a 3 V0-n 3 
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If A0 1 , A0 2 , A0 3 are the differences in <p between the two ends of the line 
segments Ax 1 , Ax 2 , Ax 3 , then 


A0‘=V0-Ax‘, A0 2 =V0-Ax 2 , A0 3 =V0 -Ax 3 


(4.27) 


If Ax 1 , Ax 2 , Ax 3 are in the directions of n 1 , n 2 and n 3 respectively, then it 
follows from equation (4.27) that 


Ap 1 „ i A c) 2 „ 2 Am 3 „ i 

— - = Vp-n , -?—=V(p- n, —A- = Vp-n- 

Ax 1 Ax 2 Ax 3 


(4.28) 


Where Ax 1 , Ax 2 and Ax 3 are the magnitudes of Ax 1 , Ax 2 and Ax 3 . Combining 
equations (4.26) and (4.28) we have 


A <p l A $? 2 


Ap 3 


Vr s= a l~T +fl 277 +a 3' 3 

Ax Ax z Ax 3 


(4.29) 


To obtain « lf « 2 a 3 , we express 
n‘=(n u n l2 n l2 ) 

n 2 = («21 «22 '*23 ) ( 4 - 3 °) 

n 3 =(n 31 n 32 n 33 ) 

Where n u , n n , n 13 are the Cartesian components of n and which can be easily 
Ax 1 Ax 1 Ax 1 

determined by — f , — f, — f. The terms, AxJ, Ax 1 ,, Ax 1 are the three 
Ax Ax Ax 

components of vector A x 1 , etc. The other values n 21 ...n 33 can be similarly 
obtained. Hence equation (4.29) can be written as 
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'n u 

n n 

n n 

r 



X' 

«21 

n 22 

n 23 


a 2 

= 

S 2 j 

_«31 

«32 

” 33 . 


a 3 




Where S {J , S 2j , S 2j are the Cartesian components of the surface vector S 
Using Cramer’s rule 


a \ 


A 

D ' 


a 2 - 


_ D 2 


D 


El 

D 


(4.31) 


Where D is the determinant of the coefficient matrix. D ] is obtained by replacing 
the first column of D by the column with the elements S }j , S 2] , and S 2j . Thus 


a { , a 2 , and a 3 are determined. 

The diffusion flux consists of two parts: normal derivative diffusion flux and 
cross derivative diffusion flux. The normal derivative flux of any dependent 
variable 0 through any cell face involves the values of 0 at cell centers whereas 
the cross-derivative diffusion flux involves the value of 0 at the edge center 
values. The normal derivative flux is treated implicitly and is coupled with the 
implicit part of the convective flux to calculates the main coefficients of the 
discretized equations while the cross-derivative diffusion flux is treated explicitly 
to avoid the possibility of producing negative coefficients in an implicit treatment. 
This term together with explicit part of convective flux is added to the source 
term. The example of the east face is taken to illustrate the diffusion model 
(Figure-4.3). Given the edge center values <p te , <p be , <p se , <p ne we can get the 


(p l? (p p ^P te ^P bs 

normal diffusion term : — , and the cross diffusion term — and 


Ax 1 


Ax 


Ax 3 


Finally, the diffusion flux is computed by 


F d - -T 

r j 


a, 


+ a. 


Ax 


+ a- 


<t>E~<! >P , „ <t>,e-<Pbe , „ <Pne-<t>se^ 

1 Ax ! 


3 Ax 3 


(4.32) 
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For the evaluation of cross-derivative flux, the following interpolation scheme is 
used to obtain the edge center values 


int * tat * trtt ' tnt 


tot 

V„ 


v TE +v P +v T +v E 


<pl 


be 


^-(pp + — <P BE + — 0 E + — (p 

v, M V V V r ‘ 


(4.33) 


(4.34) 


V lol =V BE +V P +V B +V E 


v„ 


v,„ 


vl 




Vr. 


K. 


(4.35) 


v lol =v NB +v P +v N +v E 


V, 




v F 


*«=T 


V. 


V, 


vl 


V, L. 


(4.36) 


K,,=^+V„+V J+ V £ 

Where, Vte is the volume of the cell that is located at the top-east of the cell P. 
The other edge center values of the dependent variables can be interpolated in the 
similar way. 


Sources 


The source term, which is integrated over the control volume, in the discretized 
equation can be considered as 

J v 5,dV-(s # ),V (4.37) 

Apart from the real source S^, explicitly treated parts of the convection and 
diffusion fluxes may be added to The pressure term in the momentum 
equations is also treated explicitly. Its discretization is analogous to that of the 
ordinary diffusion flux, i.e. for the i th momentum equation the pressure term is 

"J s P n f S 53 S u 

J 


(4.38) 



Where p, is the pressure at the j th face center and S tJ is the I th component of the 
surface vector for space j. 

4.7 Pressure- Velocity Coupling 


Simplified Marker and Cell (SMAC) method is used to obtain the velocity and 
pressure field that satisfies the mass and momentum equation. It is a semi-explicit 
method and it provides an efficient and easy way of pressure velocity coupling. 
The momentum equations are discritized in an explicit manner while the pressure 
gradient terms are treated implicitly, and the continuity equations are also 
enforced implicitly. Momentum and continuity equations can be expressed as 


H+I 

pV H^ + 1 { f . + f^=-'ZpT's„ 


At 

and 

£ f ; +1 =0 


(4.39) 


(4.40) 


The momentum equations are solved using the guess values of velocity and 
pressure field. The provisional velocity components U;* are calculated from the 
following equation 


Up = Up 


At 

Wr 




+ 


A t 

~pv p 


(4.41) 


Where u°p is the value of velocity at earlier iteration and S u is the pressure term. 

This provisional velocity in general will not satisfy the continuity equation. The 
continuity equation (4.40) in another form reads as follows: 

£(f;+f;) = ° 

J 



(4.42) 


j j 

Where F* is the uncorrected mass flux obtained from the provisional velocities 
and F' is the mass flux correction. To evaluate the terms on the right side of 
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equation (4.42), it is necessary that the variables (velocities and pressure) are to 
be known at the cell faces. Due to the non-staggered arrangement, if the variables 
at the cell faces are evaluated by linear interpolation between the adjacent cell 
center quantities then the pressure iteration does not converge. The solution leads 
to a Checker Board pressure field (i.e., spurious oscillation of pressure may 
occur), as shown by Rhie (1981). This problem can be circumvented by using the 
concept of Momentum Interpolation (Majumdar, 1988). The essence of the 
concept is that the velocities at the cell faces are computed by linear interpolation 
of the convective and diffusive terms but not the pressure term. Thus the method 
with the collocated variable arrangement relies indirectly on the staggering idea. 
Following this idea the interpolated velocity at the east face of the control volume 
given by is obtained in the following manner 


Where 


V P ,V E ) - — V P 

(4.43) 

p 



(4.44) 

“1- 1 

(4.45) 


And the over bar indicates a linear interpolation using equation (4.16). So the 

uncorrected mass flux for the east face using equation (4.43) becomes 

F; = pa, -S, = p(wAS, - A/Vp-S, (4.46) 


Equation (4.46) in its generalized form for any cell face is in the form 


F; = p\j- S, -ANp J -S J (4.47) 

To enforce mass conservation, the velocity and pressure corrections are 
introduced, linked by 


At 


Uj = 


Wp'i 


(4.48) 


This step simulates the effect of a staggered grid (face-center velocity and cell- 
center pressure) and reaps the benefits of faster convergence. From equation 
( 4 . 48 ) we get the mass flux correction at the face j 
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F J - P u 'j' s j = -ArVp' - S' (4.49) 

Since equation (4.49) is exactly of the same form as equation for diffusion flux 
equation (4.24), with only the variables interchanged, therefore the same method 
can be used for computing F' 

Combining equations (4.42), (4.47) and (4.49) yields 
-ArVpCS, = pV r S, -AtVp J -S J 

Or AtVp' -S, + S m = 0 (4.50) 

( s m= pVj -Sj - A tVp ] -S j) 

Which gives the pressure correction equation. After the solution of the pressure 
correction equation, the nodal velocity, mass fluxes and pressure are updated. The 
corrected pressure p n+ ' is obtained by 

// HI =//' + // (4.51) 


Instead of solving Equation (4.39) again with corrected pressure field, velocity 
correction can be calculated as follows 


u 


t __ 

p ~ 


Af 

~pV r . 




(4.52) 


The corrected velocity field u-' +1 is obtained by 


Ui n+1 = u/ + uj (4.53) 

The energy equation is discritized and solved by using this new velocity 
distribution. On integrating the conservative form of the steady state energy 
equation over a control volume, neglecting the dissipative terms, as it is 
significant for high speed flows, the temperature equations becomes 


aX(Vr-S) J =£^*-S J )r, (4.54) 

j j 

which is solved by SOR method to get the temperature distribution at the current 
time step. 
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The first terms of equation (4.50) and equation (4.54) are expanded following the 
philosophy that was followed during the discretization of equation (4.24). 
Subsequent regrouping of the normal derivative terms yields the standard Poisson 
Equation. The Poisson equation can be solved by a multitude of methods. We 
have used the Gauss-Siedel iterative method, which is quite reliable and robust. 
The flow chart with the main steps of the algorithm is given in Figure 4.4. 


4.8 Solution Algorithm 

The velocity pressure and temperature fields are calculated by the following the 
procedure of Eswaran and Prakash (1998). The procedure can be summarized in 
the following way: 

1. Make initial guess of the pressure, velocity and temperature fields. Use 
equations (4.41) and (4.44) to calculate cell center u P and v P . Use linear 
interpolation equation (4.16) to obtain the face center quantities v,. 

2. Compute the mass flux through cell face j using equation (4.47). 

3. Use equation (4.48) to compute the flux correction at the cell face j using 
equation (4.49), which is computed using the formulation for diffusion fluxes 
[ Equations (4.24) to (4.36) ] For calculating by Fj , the variable <j) is to be 

replaced by p . 

4. Compute the residue, 91 for each cell 

91 =-£ F,‘-£ f,' 

J J 

5 Calculate the cell-center pressure correction (using Gauss-Siedel method) 
from the relation 
91 

p » =Pp +C 0 ^a? 

where co is the relaxation factor and ap stands for the diagonal coefficients of 
equation (4.50), which is calculated from 


a x 

gi 

. a 2 1 

a 2 

a, 

4- J - 

a 3 
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Ax 2 
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Ax 2 

, Axi 

Ax 3 
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where a l ,a 2 etc are the same as that are used in equation (4.24) to (4.36). 

6 If 9t mf > e go to step 3. 

7 Store the updated mass flux through cell faces using equation (4.47). 

8 Store the updated pressure at the center of the cells 
Pp = Pp + Pp 

9 Store the cell-center corrected velocity 

At 


Up — 


n+1 


pV, 


YjP'j S ‘J 


p J 


u, ' = u, + u, 

10 Calculate the residual for each cell 

$ = «£(vr.s),=2>s ; )7-, 

J J 

1 1 Calculate the cell-center temperature using Gauss-Seidel method to obtain 

9t 


Tp = Tp + S2 


b p a 


Where Q is the relaxation factor and bp stand for the diagonal coefficients of 
equation (4.54) i.e., calculated from the relation 

bp = - cip 

12 If > e go to step 10. 

13 Go to step 1 and repeat the process until the steady state is reached. It can be 
shown that solving equation (4.42) is same as solving 

V-u,* 


V 2 p' =-p- 

and finding the corrections 

Vf • 
u = Vp 

J P 


At 


such that 


n+l * , / 

u, = u. + u, 


p " +1 = p ” + p ' 
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Thus the procedure above can be called solution of Pressure Correction 
Linked Equations. 


4.9 Numerical Stability Considerations 


Semi-Explicit method used for present computations relies on explicit differences 
and impose certain restrictions on time step. For a given mesh the choice of time 
step is determined through stability analysis, which has to take care of two 
conditions. First, fluid should not be allowed to cross more than one cell in any 
one-time step. This restrictions is derived from the Courant-Friedrichs-Lewy 
(CFL) conditions given by 


<5ti < min 


Sx. Sy 5z 

U ’ V ’|w| 


(4.55) 


Where the minimum is with respect to every cell in the domain. Typically S ti is 
chosen equal to one-third to two-third of the minimum cell transient time. 
Secondly, when a non-zero value of kinematic viscosity is used momentum must 
not diffuse more than approximately one cell in one time step. A linear stability 
analysis shows that the restrictions on grid Fourier number will yield 


1 (£x) 2 (<fy) 2 (&) 

2(<5x) 2 +(<5y) 2 +(&) 2 


(4.56) 


Finally, the minimum of the two above time increments is chosen for 
computation. 
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Figure 4.1 Schematic of grid mesh and computational domain 




Figure 4.3 Face representation to illustrate the diffusion model 




Figure 4.4 


Flowchart of iterative solution scheme 










Chapter 5 


Results and Discussion 


5.1 Introduction 

A three dimensional numerical solution of flow and heat transfer for a plane 
channel with built in circular tube and delta winglet pair in common-flow-up 
configuration is presented in this section. A 74x48x23 grid-mesh is used in the 
present computation. The divergence-free criterion is satisfied using an upper 
bound of lO^ 4 . Computations have been carried out for moderately high Reynolds 
numbers in the range of 1000 and 1350. Air has been taken as the working fluid; 
hence the Prandtl number of the present study is 0.7. The Nusselt number is 
calculated on the basis of bulk mean temperature. D/B has been kept constant 
equal to 0.36 for all computations. 

5.2 Flow and Heat Transfer Characteristics in a Rectangular Channel with a 
built-in Circular Tube 

A detailed analysis of flow and heat transfer characteristics in a rectangular 
channel with a built-in circular tube is presented herein. The results presented 
here are for a Reynolds number of 1000. 

5.2.1 Flow Characteristics 

Figure 5.1 (a) shows limiting streamlines near the bottom wall for the 
instantaneous flow field in a channel with a built-in circular tube. A saddle point 
of separation and a horseshoe vortex system is observed in the figure. The 
incoming flow reaches a stagnation or saddle point of separation (marked by A) 
and goes around the body. The nodal point of attachment (marked as C) and the 
separation lines are also visible. On each side of the tube, one finds a region of 


41 



converging streamlines (marked as G). These are the traces of horseshoe vortices. 
The symbol E marks the footprints of shed vortices. Figure 5.1 (b) shows the 
limiting streamlines close to the bottom for the time-averaged flow field. The 
traces of symmetric vortices are seen in the figure. 

Figures 5.2 (a) through 5.2-(f) show vortex shedding past the tube at the 
horizontal midplane of the channel. The sequence of events in the figure reveals 
that the vortex from the bottom side of the tube is shed into the stream while 
another vortex grows at that same location to the limit of onset of shedding. 
Likewise, Figures 5.3 (a) through 5.3-(f) show the shedding of the vortex from the 
topside of the tube. Figure 5.4 shows the streamlines for the time-averaged flow 
on the horizontal midplane. 

5.2.2 Heat Transfer Characteristics 

Figure 5.5 represents the Nusselt number distribution on the bottom plate of the 
channel with a built-in circular tube for the Reynolds number of 1000. The 
following regimes may be distinguished. At the leading edge of the bottom plane, 
Nusselt number starts with a high value and decreases gradually. At the leading 
edge of the plate, the cooler fluid comes in contact with the hot fin for the first 
time hence the heat transfer is very high. The gradual decrease in the Nusselt 
number is attributed due to the development of thermal boundary layer on the 
channel walls. In front of the tube, a zone of high Nusselt Number is observed. 
This increase in Nusselt number results from the formation of the horseshoe 
vortex system that consists of two counter-rotating longitudinal vortices. As the 
fluid approaches the stagnation line of the circular tube, it slows down and its 
pressure increases. The smaller velocity within the boundary layer in the vicinity 
of the bottom plate, which supports the circular tube, leads to a smaller pressure 
increase. Thus, the induced pressure gradient on the stagnation line causes the 
flow towards the bottom wall that interacts with the main stream. The fluid rolls 
up forming vortices, which finally wraps around the front half of the tube and 
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extends to the rear of the tube (Goldstein and Kami (1984)). The horseshoe 
vortices generated in front of the tube bring about better mixing and the heat 
transfer in this region is enhanced significantly. The Nusselt number is low in the 
wake region .The poor heat transfer in this region is attributed to the separated 
dead water zone with fluid recalculating at a low velocity. 

5.3 Effect of Winglet-Pair on Flow and Heat Transfer 

Detailed analysis of flow structure and heat transfer characteristics in a channel 
with a built-in tube along with winglet pair placed in common-flow-up 
configuration is presented herein. The results presented correspond to a Reynolds 
number of 1000. Different values of transverse distance between the trailing edges 
of the winglets have been considered. However, generally b/D =1.08 for delta 
winglet and b/D =1.67 for rectangular winglet, if it is not mentioned otherwise in 
the text. 

5.3.1 Flow Characteristics 

Figure 5.6 shows the limiting streamlines for the velocity field on a horizontal 
plane close to the bottom wall with winglets placed in common-flow-up 
configuration behind the tube. If one considers the flow passage between one of 
the winglets and the curved surface of the tube, the geometrical shape resembles a 
nozzle like configuration. The fluid layers do not face adverse pressure gradient in 
the aft of the circular tube. Instead, the fluid accelerates in the constricted passage 
and the wake region becomes quite narrow. The winglet pair causes transverse 
motion, which in turn delays the separation. Separation is not observed till a 
separation angle of 170°. A very narrow wake is formed indicating F as the wake 
stagnation point. 

Figure 5.7 shows streamline plot on the horizontal mid-plane with winglets in 
common-flow-up configuration. Symmetry can be clearly observed even in an 
instantaneous flow field in this case as vortex shedding is suppressed. 
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Figure 5.8 shows streamlines on the vertical mid-plane with winglets placed 
behind tube in common-flow-up configuration. The formation of horseshoe 
vortices in front of the tube is evident from this figure. Placing of winglets stops 
the vortex shedding and creates a strong swirling motion behind the tube. The 
swirl is primarily characterized by a longitudinal vortex system. The projection of 
this complex structure is observed on the vertical plane behind the tube. 

Figures 5.9 (a), (b) and (c) show cross-stream velocity vectors at various axial 
locations in the channel with the circular tube and the winglet pair as the built-in 
obstacles. Schematic location of the winglet pair has already been shown in 
Figure 1.2. The axial location of Fig. 5.9(a) is near the trailing edge of the winglet 
pair. The velocity vectors depict the secondary flow pattern on various cross- 
stream planes in the channel. The results correspond to a Reynolds number of 

1000 and an angle of attack of 138.61 0 . The pressure difference between the front 
surface, facing the flow and the suction surface of the winglet generates swirling 
motion. The swirling motion is described by the longitudinal vortices. The figure 
clearly reveals that the crossflow pattern depicts a common-flow-up situation. 
Mendez et al. (1998) have explained the kinematics of flow and the mechanism of 
transport enhancement due to longitudinal vortices. The strength of the vortex 
motion decreases with the downstream distance from the winglet. The decrease in 
vortex strength can be attributed to viscous diffusion. 

Figures 5.10 (a), (b) and (c) show cross-stream velocity vectors at same axial 
locations in the channel but with rectangular winglet of same aspect ratio (double 
the area of the delta winglet), instead. At the same axial location, the vortex 
strength described by the rectangular winglet is higher as compared to the delta 
winglet. The surface area of the rectangular winglet is more than that of the delta 
winglet. 

Figures 5.11(a), (b), (c) and Figures 5.12(a), (b), (c) represent the streamlines on 
the cross-stream plane at the non-dimensional axial locations of 10.17, 13.75 and 
17.32 from the leading edge of the channel for the delta and rectangular winglets 
respectively. The streamlines confirm existence of a strong vortical motion that 


44 



was explained through Figures 5.9 and 5.10. The centers of the main vortices of 
the longitudinal vortex system are located at distance of ± 0.3B from the vertical 
midplane. Induced vortices become prominent as the pair of longitudinal vortices 
travel in the downstream. 

The lateral distance between the winglet pair has been varied to observed its 
effect on the vortical motion. The lateral distance between the trailing edges is 
given by b. For the cases, described in Figures 5.9 and 5.10, the parameter b ID is 
taken as 1.67. For the delta winglet pairs and b/D of 1.08, Figs 5.13 (a), (b) and 
(c) have described the vortical motion at three different cross-stream planes. One 
can see that the main vortices of the longitudinal vortex system have come closer 
to each other due to reduced b/D. 

Figures 5.14, 5.15 and 5.16 show three-dimensional particle paths for the delta 
and rectangular winglets placed in common-flow-up configuration. Complex 
vortex structures generated due to the configurations of interest are clearly visible 
in these figures. All the figures show generation of horseshoe vortices in front of 
the tube. The generation longitudinal vortices due to the delta winglet pair with 
b/D = 1.67 has been shown in Figure 5.14. Likewise, Fig. 5.15 shows the 
generation of longitudinal vortices due to the rectangular winglet pair. Figure 5.16 
displays the vortical structure due to delta winglet pair with b/D = 1.08. 

5.3.2 Heat Transfer Characteristics 

Nusselt number distribution at bottom plate of channel with winglets placed in 
common-flow-up configuration is shown in Figure 5.17. From entrance of the 
channel, till the axial location of the leading edge of winglet pair, the variation of 
Nu is same as Figure 5.5. Beyond this location, the Nussselt number starts 
increasing. A sharp rise in Nusselt number is observed near the location of the 
trailing edge of the winglet. This enhancement is attributed to the formation of a 
complex vortex system that consists of two counter-rotating longitudinal vortices. 
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The formation of longitudinal vortices brings about a better mixing and the 
Nusselt number becomes high in the zone influenced by the swirl. 

Figure 5.18 compares the span-averaged Nusselt number distribution on the 
bottom wall in a channel for the cases of built-in circular tube, without and with 
winglet pair. The Reynolds number of interest is 1000. The span-averaged Nusselt 
number reaches the maximum value near the axial location of the leading edge of 
the circular tube. This is primarily governed by the formation of horseshoe 
vortices at the tube and bottom wall junction as explained earlier. The peak value 
of the span-averaged Nusselt number is about 11.2. Up to the location of the 
leading edge of the winglets, the distribution of span-averaged Nusselt number is 
almost same for the cases with and without winglets. However, beyond that 
location, Nusselt number starts increasing for the case of built-in winglet pair. The 
maximum enhancement in heat transfer (95 percent) takes place at the immediate 
upstream of the trailing edge of the winglet pair (X = 9.8). At the exit of the 
channel, the span-averaged value of Nusselt number is 5.8. The span-averaged 
Nusselt number value at the exit (X = 20) of a plane channel is 4.0 while for the 
case of a channel with a built-in circular tube, it is 3.7. 

Figure 5.19 compares the distribution of span-averaged Nusselt number on the 
bottom wall of channel for two Reynolds numbers, namely Re = 1000 and Re = 
1350. As expected, the Nusselt number distribution is higher for the higher 
Reynolds number. 

Figure 5.20 compares the distribution of span-averaged Nusselt number on the 
bottom wall of the channel with a built-in circular tube for cases of delta and the 
rectangular winglet pair as the second obstacle. The Nusselt number distribution 
is higher for the rectangular winglet pair as compared to the delta winglet pair. 
The surface area of the rectangular winglet pair being more (due to same aspect 
ratio of two different types of winglets), the strength of the vortical motion in the 
case of rectangular winglet pair is high. 

Figure 5.21 shows the distribution of span-averaged Nusselt number on the 
bottom wall of the channel with built-in circular tube for the cases of delta 
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Abstract 


Numerical investigation of flow and heat transfer in a rectangular channel with 
built-in circular tube and different types of winglets placed in common-flow-up 
configuration has been carried out. The Reynolds numbers are moderately high, in 
the range of 1000. Existing air-cooled condensers in geothermal power plants use 
fin-tube heat exchangers with circular tubes. The size of the heat exchangers is 
huge and often the cost of the condensers is more than one-third of the plant cost. 
The size of the condensers can be reduced through enhancement of heat transfer 
from fin surfaces. The enhancement strategy involves introduction of strong 
swirling motion in the flow field. The swirl can be generated by the longitudinal 
vortices. In this study, the longitudinal vortices are created by the winglet type 
vortex generators, which are mounted behind the tubes. An element of a heat 
exchanger has been considered for detailed study of the flow structure and heat 
transfer analysis. 

In the present study, a novel technique (Torii et. al., 2000) has been utilized for 
the enhancement of transport coefficients. The winglets are placed with a 
heretofore-unused orientation for the purpose of augmentation of heat transfer. 
This orientation is called as common-flow-up configuration. The proposed 
configuration causes significant separation delay, reduces form drag, and removes 
the zone of poor heat transfer from the near-wake of the tubes. The analyses of 
flow and heat transfer in the proposed configuration have been accomplished 
through a numerical solution of complete Navier-Stokes and energy equations. 



winglets with different transverse distances between their trailing edges. The 
span-averaged Nusselt number distribution is almost same for the two situations, 
mentioned above. For the case of lateral distance b=1.08D, the span-averaged 
Nusselt number increases in the wake region of the winglets, whereas the span- 
averaged Nusselt number is slightly less near the leading edge as compared to the 
case of b=1.67D. 

5.4 Spatial Grid Independence Model Validation 

Spatial Grid independence test has been carried out using two grid-meshes of 
sizes 74 x48x23 and 60x48x20. Figure 5.22 shows span-wise average Nusselt 
number distribution for two grid-meshes. They compare favorably indicating the 
grid-independent situation for the grid size of 74 x48x23. 

The code validation was accomplished through comparison with the experimental 
results of Torii (2001). Figure 5.23 compares the present computation with the 
experimental results of Torii (2001). The enhancement is expressed in terms of 

the ratio of the ratio of Nu s and Nu s0 . Here Nu s is the span-averaged Nusselt 
number at any axial location for the case of a channel with built-in circular tube 
and winglet pair and Nu s0 is the span-averaged Nusselt number at any axial 
location for the case of a channel with a built-in circular tube. The axial location 
of the center of the tube has been taken as the reference coordinate point in the 
figure. The Reynolds number of interest is 1320 and the experimental technique 
of heat transfer measurement is based upon the experimental study of Mochizuki 
et al (1988). The computation and experimental results compare well in general, 
and in the downstream of the trailing edge in particular. 
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Figure 5.9 Cross-stream velocity vectors at different x-locations behind the delta 
winglet-pair (b/D=1 .67) 
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Figure 5.10 Cross-stream velocity vectors at different x-locations behind the rectangular 
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Figure 5.12 Streamlines on the cross-stream planes at different x-locations behind the 
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Figure S.13 Streamlines on the cross-stream planes at different x-locations behind the 
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Figure 5.14 Three dimensional particle paths behind the delta winglet pair in common- 
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Figure 5.16 Three dimensional particle paths behind the delta winglet-pair in common- 
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Figure 5.18 Nu s for a plane channel, channel with built-in circular-tube and channel with 
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Chapter 6 


Conclusion and Scope for Future Work 

6.1 Concluding Remarks 

A computational study has been accomplished to determine the flow structure and 
heat transfer in a rectangular channel with a built-in circular tube and different 
shaped vortex generators. The horseshoe vortices are generated at the junction of 
the forward stagnation line and the bottom plate. The horseshoe vortices travel in 
the downstream inducing a flow structure, which is primarily a longitudinal 
vortex system. However, this natural horseshoe vortex system does not persist in 
the downstream. Stronger longitudinal vortices are created by placing the winglet 
type vortex generators behind the circular tube. The present analysis reveals that 
the combination of circular tube and the winglet pair with common-flow-up 
configuration significantly improves the heat transfer performance. The nozzle 
like flow passages created by the winglet pair and the aft region of the circular 
tube promote accelerating flow and there by remove the zone of poor heat transfer 
from the near wake. At the end of a channel (X=20) the span-averaged value of 
Nusselt number (based on channel height) is 5.8 for a Reynolds number of 1000. 
It may be mentioned that span-averaged Nusselt number (based on channel 
height) near the exit of a plane channel is 4.0 for the same Reynolds number. 
Enhancement of heat transfer is remarkable even at far downstream. The average 
Nusselt number of the entire bottom plate for the case of a channel with built in 
circular tube is 5.37 for a Reynolds number of 1000. The average Nusselt number 
of the plate improves to a value of 7.27 for the same Reynolds number when the 
winglet pair is added following the common-flow-up configuration. 
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6.2 Scope for Future Work 


The results of this work reveal that the vortex generators with a common-flow-up 
configuration as a useful device for improving heat transfer in the air-cooled 
condensers for the geothermal power plants and other heat exchangers. The oval 
tubes can replace the circular tubes. This will reduce the form drag. 

Here the computations have been done assuming flow regime to be laminar. 
Though turbulent flow is not frequently encountered in fin-tube heat exchangers, 
some special applications involving very high velocities, the flow regime can 
become turbulent. Therefore present study can be extended for the turbulent 
flows. Using an appropriate turbulence model, the performance of the proposed 
design can be computed for higher Reynolds numbers. 


73 



References 


1. Biswas, G. and Chattopadhaya, H. (1992), Heat Transfer in a Channel with Built- 
in Wing-Type Vortex Generators, bit. J. Heat Mass transfer , Vol. 35, pp. 803- 
814. 

2. Biswas, G., and Mitra, N. K. and Fiebig, M. (1994), Heat Transfer Enhancement 
in Fin-Tube Heat Exchangers by Winglet-Type Vortex Generators, bit. J. Heat 
Mass Transfer, Vol. 37, pp. 283-291. 

3. Biswas, G., Deb, P., and Biswas, S., (1994), Generation of Longitudinal 
Streamwise Vortices - A Device for Improving Heat Exchanger Design, Journal 
of Heat Transfer (ASME), Vol. 116, pp. 588-597. 

4. Biswas, G., Torii, K., Fujii, D. and Nishino, K. (1996), Numerical and 
Experimental Determination of Flow Structure and Heat Transfer Effects of 
Longitudinal Vortices in a Channel Flow, Int. J. Heat Mass Transfer, Vol. 39, No. 
16, pp. 3441-3451. 

5. Cantwell, B. and Coles, D., (1983), An Experimental Study of Entrainment and 
Transport in the Turbulent Near Wake Circular Cylinder, J. Fluid Meek, Vol. 
136, pp. 321-374. 

6. Eibeck, P.A. and Eaton, J.K. (1987), Heat Transfer of a Longitudinal Vortex 
Embedded in a Turbulent Boundary Layer, Journal of Heat Transfer (ASME), 
Vol. 6, pp. 16-24. 


74 



7. Eswaran, V. and Prakash, S. (1998), A Finite Volume Method for Navier-Stokes 
Equations, Proceedings of the third Asian CFD Conference, Bangalore, Vol. 1, 
pp. 251-269. 

8. Feigenbaum, M. J., (1980), The Onset Spectrum of Turbulence, Phys. Lett., part- 
A, Vol. 74, pp. 375-378. 

9. Ferziger, J. H., and Peric, M. (1999), Computational Methods for Fluid 
Dynamics, Springer Verlag, Berlin, New York 

10. Fiebig, M., Kallweit, P. and Mitra, N. K. (1986), Wing Type Vortex Generators 
for Heat Transfer Enhancement, Proc. 8 th 1HTC, Vol. 6, pp. 2909-2913. 

11. Fiebig, M., Brockmeier, U., Mitra, N. K. and Guntermann, T. (1989), structure of 
Velocity and Temperature Fields in Laminar Channel Flows with Longitudinal 
Vortex Generators, Numerical Heat Transfer, Part A, Vol. 15, pp. 281-302 

12. Fiebig, M., Dong, Y. and Mitra N.K. (1990), Warmeuberganangserhohung und 
Widerstandsverringerung durch Laenhgswirbelerzeuger in Rippenrohren, Waeme 
and Stoffubertragung, Vol. 25, pp. 34-44. 

13. Fiebig, M., Kallweit, P., Mitra, N. K. and Tiggelbeck, S. (1991), Heat Transfer 
Enhancement and Drag by Longitudinal Vortex Generators in Channel Flow, 
Experimental Thermal and Fluid Science, Vol. 4, pp. 103-114. 

14. Garg, V. K. and Maji, P. K. (1987), Flow through a Converging-Diverging Tube 
with Constant Wall Enthalpy, Numerical Heat Transfer , Vol. 12, pp. 285-305. 


75 



15. Gentry, M. C., Jacobi, A. M. (1997) , Heat Transfer Enhancement by Delta-Wing 
Vortex Generators on a Flat Plate: Vortex Interactions with Boundary Layer, 
Experimental Thermal and Fluid Science, Vol. 14, pp. 231-242 

16. Goldstein, R. J. and Kami, J.(1984), The Effect of a Wall Boundary Layer on a 
Local Mass Transfer from a Cylinder in Crossflow, Journal of Heat Transfer 
(ASME), Vol. 106, pp.260-267 

17. Guzm&n, A. M. and Amon, C. H., (1996), Dynamical Flow Characterization of 
Transitional and Chaotic regimes in Converging-Diverging Channels, J. Fluid 
Mech., Vol. 321, pp. 25-57. 

18. Harlow, F. H., and Welch, J. E., (1965), Numerical Calculation of Time- 
dependent Viscous Incompressible Flow of Fluid with Free Surface, The Phys. 
Fluids , Vol. 8, pp. 2182-2188. 

19. Hirt, C. W. and Cook, J. L. (1972), Calculating Three-Dimensional Flows Around 
Structures and Over Rough Terrain, J. Comp. Phys., Vol. 10, pp. 324-340. 

20. Khosla, P. K., and Rubin, S. G., (1974), A diagonally dominant second-order 
accurate implicit scheme, Comp. & Fluids, Vol. 2, pp. 207-209. 

21. Kordulla, W., and Vinokur, M., (1983), Efficient Computation of Volume in Flow 
Predictions, AIAA J., Vol. 21, No. 6, pp. 917-918. 

22. Kobayashi, M.H. and Pereira, C.F. (1991), Calculation of Incompressible Laminar 
Flows on a Non-Staggered, Non-Orthogonal Grid, Numerical Heat Transfer, Part 
B, Vol. 19, pp. 243-262. 


76 



4.5 Surfaces areas and Volumes 23 

4.6 Discretization Procedure 23 

4.6.1 Discretization of Continuity Equation 24 

4.6.2 Discretization of General Equation 24 

4.7 Pressure- Velocity Coupling 31 

4.8 Solution Algorithm 34 

4.9 Numerical Stability Considerations 36 

5 Results and Discussion 

5.1 Introduction 41 

5.2 Flow and Heat Transfer Characteristics in a Rectangular Channel 41 

with a built-in Circular Tube 

5.2.1 Flow Characteristics 41 

5.2.1 Heat Transfer Characteristics 42 

5.3 Effect of Winglet-Pair on Flow and Heat Transfer 43 

5.3.1 Flow Characteristics 43 

5.3.2 Heat Transfer Characteristics 45 

5.4 Spatial Grid Independence Model Validation 47 

6 Conclusion and Scope for Future Work 

6.1 Concluding Remarks 72 

6.2 Scope for future work 73 

7 References 74 



23. Majumdar, S., Rodi, W., and Zhu, J., 1992, Three-Dimensional Finite-Volume 
Method for Incompressible flows with Complex Boundaries, J. of Fluids Engg. 
(ASME), Vol. 114, pp. 496-503. 


24. Majumdai, S. (1988), Role of Undei relaxation in Momentum Interpolation for 
calculation of Flow with Non-Staggered Grids, Numerical Heat Transfer, Vol. 13, 
pp. 125-132. 

25. Manneville, P. and Pomeau, Y. (1980), Different ways to Turbulent in Dissipative 
Dynamical Systems, Physica D, Vol. 1, pp. 219-226. 

26. Mendez, R. R., Sen, M., Yang, K. T. and McClain, R. L. (1998), Enhancement of 
Heat transfer in an Inviscid-Flow Thermal Boundary layer Due to a Rankine 
Vortex, Int. J. Heat Mass Transfer, Vol. 41, pp. 3829-3840. 

27. Mochizuki, S., Yogi Y. and Yang W. J. (1988), Advances in Single-Blow Method 
for Performance Evaluation of Heat Transfer Surfaces, Proc. of 2 nd Int. Symp. On 
Heat Transfer, pp. 284-291, Beijing, China. 

28. Mukhopadhyay, A., Sundarajan, T. and Biswas, G. (1993), An Explicit Transient 
Algorithm for Predicting Incompressible Viscous Flows in Arbitrary Geometry, 
Int. J. of Heat Transfer, Vol. 114, pp. 280-284. 

29. Nichols, B. D. and Flirt, C. W., (1971), Improved Free Surface Boundary 
Conditions for Numerical Incompressible Calculations, J. Compt. Phys., Vol. 8, 
pp. 434-448. 

30. Orlanski, I. (1976), A Simple Boundary Condition for Unbounded Flows, J. 
Comp. Phys, Vol. 21, pp. 251-269. 


77 



31. Patankai, S. V. and Spalding, D. B. (1972), A Calculation Procedure for Heat 
Mass and Momentum Transfer in three-Dimensional Parabolic Flows, Int. J. 
Heat Mass Transfer, Vol. 15, pp. 1787-1806. 

32. Patankai, S. V. (1981), A Calculation Procedure for Heat and Momentum 
Transfer in Three-Dimensional Parabolic Flows, Int. J. Heat Mass Transfer, Vol. 
15, pp. 1787-1806. 

33. Pauley, W. R. and Eaton, J. K. (1988), Experimental Study of the Development of 
Longitudinal Vortex pairs Embedded in a Turbulent Boundary Layer, AIM J., 
Vol. 26, pp. 816-823. 

34. Peric, M. (1985), A Finite Volume Method for the Prediction Of Three- 
Dimensional Fluid Flow in Complex Ducts, Ph. D. Thesis, University of London. 

35. Peric, M., Kessler, R. and Scheuerer, G. (1988), Comparison of Finite-Volume 
Numerical Methods with a Staggered and Collocated Grids, Computers & Fluids, 
Vol. 16, pp. 389-403. 

36. Pulliam, T. H. and Vastano, J. A., (1993), Transition to chaos in an Open 
Unforced 2D Flow, J. Compt. Phys., Vol. 96, pp. 15-53. 

37. Robichaux, J., Tafti, D. K. and Vanka, S. P., (1992), Large Eddy Simulations of 
Turbulence on the CM-2, Numerical Heat Transfer, Part-B, Vol. 21, pp. 367-388. 

38. Ruelle, D. and Takens, F., (1971), On the Nature of turbulence, Commun. Math. 
Phys., Vol. 20, pp. 167-192. 

39. Rhie, C. M. and Chow, B. L. (1983), Numerical Study of Flow Past an Aerofoil 
with Trailing Edge Separation, AIAA J., Vol. 21, pp. 1525-1532. 


78 



40. Rhie, C. M. (1981), A Numerical Study of the Flow Past an Isolated Airflow with 
Separation, Ph.D. Thesis, Dept, of Mechanical and Industrial Engg., University of 
Illinios at Urbana Champaign, 1981. 

41. Sanchez, M., Mitra, N. K. and Fiebig, M. (1989), Numerical Investigations of 
Three-Dimensional Laminar Flows in a Channel with Built-In Circular Cylinder 
Wing-Type Vortex Generators. Proc. Eighth GAMM-Conference on numerical 
methods, Vieweg Verlag, pp. 484-492 


42. Tiggelbeck, S., Mitra, N. K., and Fiebig, M., (1994), Comparison of Wing Type 
Vortex Generators for Heat Transfer Enhancement in Channel Flows, Journal of 
Heat Transfer (ASME), Vol. 116, pp. 880-885. 

43. Torii K., Nishino K. M., Kwak K. M., and Kawai R. (2000), YNU Team Report: 
Research Results and Implementation Plan, NEDO VORTEX Team Meeting. 

44. Torii K., (2001) Personal Communications. 

45. Turk, A. J. and Junklan, G. H. (1986), Heat Transfer Enhancement Downstream 
of Vortex Generators on a Flate Plate, Proc. 8 ,h IHTC, Vol 6 , pp. 2903-2908, San 
Francisco 

46. Tsai, S. F., Sheu, T. W. FI. and Lee, S. M. (1999), Heat transfer in a conjugate 
heat exchanger with a wavy fin surface, Int. J. Heat Mass Transfer, Vol. 42, pp. 
1735-1745. 

47. Valencia, A., Fiebig, M., and Mitra, N. K., (1996) Heat Transfer Enhancement by 
Longitudinal Vortices in Fin-Tube Heat Exchanger Element with Flat Tubes, 
Journal of Heat Transfer (ASME), Vol. 1 18, pp 209-211. 


79 



48. Van Doormal, J. P., and Raithby, G. D., (1984), Enhancement of SIMPLE method 
of Predicting Incompressible Fluid Flows, Numerical Heat Transfer, Vol. 7, pp 
147-163. 

49. Veima, A. K. and Eswaran, V.(1996), Overlapping Control Volume Approach for 
Convection-Diffusion Equation, Int. J. Numerical Methods Fluids, Vol. 23, pp. 
865-882. 


50. Verma A. K., and Eswaran, V.(1997), A Bounded Convection Scheme for the 
Overlapping Control Volume Approach, Int. J. Numerical Methodes Fluids 

51. Viltori, G. and Blondeaux, P. (1993), Quasiperiodicity and Phase Locking Route 
to chaos in the Two-Dimensional Oscillatory Flow Around a Circular Cylinder, 
Phys. Fluids, Part-A, Vol. 5, pp. 1866-1868. 

52. Westphal, R. V., and Mehta, R. D., 1987, Interaction of Oscillating Vortex with a 
Turbulent Boundary layer, AIAA 19 th Fluid Dynamics and Laser Conference, 
Hawaii, AIAA paper no. 87-1246. 

53. Yanagihara, J. L. and Torii, K., (1990), Heat Transfer characteristics of Laminar 
Boundary Layer in the Presence of Vortex Generators, Proc9 th IHTC, Vol. 6, 
pp.323-328. 


80 




List of Figures 


Figure 1 . 1 Schematic diagram of core region of a fin-tube heat exchanger 4 

Figure 1.2 Delta-winglet type vortex generators on flat surface 4 

Figure 1.3 Heat exchanger module considered in present investigation 5 

Figure 2.1 Element of existing fin-tube heat exchanger 14 

Figure 4.1 Schematic of the grid-mesh and the computational domain 37 

Figure 4.2 Three-dimensional finite volume cell 38 

Figure 4.3 Face representation to illustrate diffusion model 39 

Figure 4.4 Flowchart of iterative solution scheme 40 

Figure 5.1(a) Limiting streamlines based on instantaneous flow on a plane 

close to the bottom plate for a channel with built-in circular tube 48 
(b) Limiting streamlines based on time-averaged flow on a plane 

close to the bottom plate for a channel with built-in circular tube 49 
Figure 5.2 Vortex shedding at the mid-plane in the wake of a circular tube 

(shedding starts from the bottom side of the tube). 50 

Figure 5.3 Vortex shedding at the mid-plane in the wake of a circular tube 

(shedding starts from the top side of the tube). 51 

Figure5.4 Streamline-plot for time-averaged flow on a horizontal mid-plane 

of the channel with a built-in circular tube 52 



